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The effect of the lateral interactions on the critical behavior of long 
straight rigid rods on two-dimensional lattices 

P. LongonePD. H. LinaresP A. J. Ramirez-Pastor^* 

Using Monte Carlo simulations and finite-size scaling analysis, the critical behavior of 
attractive rigid rods of length k (fc-mers) on square lattices has been studied. An ordered 
state, with the majority of fc-mers being horizontally or vertically aligned, was found. This 
ordered phase is separated from the disordered state by a continuous transition occurring at 
a critical density 6 C , which increases linearly with the magnitude of the lateral interactions. 



I. Introduction 

The study of systems of hard non-spherical col- 
loidal particles has, for many years, been attract- 
ing a great deal of interest and the activity in this 
field is still growing [1-14]. An early seminal con- 
tribution to this subject was made by Onsager [1] 
with his paper on the isotropic-nematic (I-N) phase 
transition in liquid crystals. The Onsager's theory 
predicted that very long and thin rods interacting 
with only excluded volume interaction can lead to 
long-range orientational (nematic) order. Thus, at 
low densities, the molecules are typically far from 
each other and the resulting state is an isotropic 
gas. However, at large densities, it is more favor- 
able for the molecules to align spontaneously (there 
are many more ways of placing nearly aligned rods 
than randomly oriented ones), and a nematic phase 
is present at equilibrium. 

Interestingly, a number of papers have appeared 
recently, in which the I-N transition was studied 
in two dimensions [10-14]. In Ref. [TU], the au- 
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thors gathered strong numerical evidence to sug- 
gest that a system of square geometry, with two 
allowed orientations, shows nematic order at inter- 
mediate densities for k > 7 and provided a qualita- 
tive description of a second phase transition (from 
a nematic order to a non-nematic state) occurring 
at a density close to 1. However, the authors were 
not able to determine the critical quantities (crit- 
ical point and critical exponents) characterizing 
the I-N phase transition occurring in the system. 
This problem was resolved in Refs. [TTJ [T2] , where 
an accurate determination of the critical expo- 
nents, along with the behavior of Binder cumulants, 
showed that the transition from the low-density dis- 
ordered phase to the intermediate-density ordered 
phase belongs to the 2D Ising universality class for 
square lattices and the three-state Potts univer- 
sality class for honeycomb and triangular lattices. 
Later, the I-N phase transition was analyzed by 
combining Monte Carlo (MC) simulations and the- 
oretical analysis [Ml [14]. The study in Refs. [LB) [14] 
allowed (1) to obtain 8 C as a function of k for 
square, triangular and honeycomb lattices, being 
9 c (k) oc k" 1 (this dependence was already noted 
in Ref. [10] ) : and (2) to determine the minimum 
value of k (k m i n ), which allows the formation of a 
nematic phase on triangular {kmin = 7) and hon- 
eycomb (k m i n — 11) lattices. 
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In a recent paper, Fischer and Vink [15) indi- 
cated that the transition studied in Refs. [10-14] 
corresponds to a liquid-gas transition, rather than 
I-N. This interpretation is consistent with the 2D- 
Ising critical behavior observed for monodisperse 
rigid rods on square lattices |llj . This point will 
be discussed in more detail in Sec. III. 

In contrast to the systems studied in Refs. [10- 
14], many rod-like biological polymers are formed 
by monomers reversibly self-assembling into chains 
of arbitrary length. Consequently, these systems 
exhibit a broad equilibrium distribution of fila- 
ment lengths. A model of self- assembled rigid 
rods has been recently considered by Tavares et 
al. [16] . The authors focused on a system composed 
of monomers with two attractive (sticky) poles that 
polymerize reversibly into polydisperse chains and, 
at the same time, undergo a continuous I-N phase 
transition. The obtained results revealed that ne- 
matic ordering enhances bonding. In addition, the 
average rod length was described quantitatively in 
both phases, while the location of the ordering tran- 
sition, which was found to be continuous, was pre- 
dicted scmiquantitatively by the theory. 

Beyond the differences between lattice geometry 
and the characteristics of the rods (self-assembled 
or not), one fundamental feature is preserved in 
all the studies mentioned above. This is the as- 
sumption that only excluded volume interactions 
between the rods are considered (except in Ref. [16] , 
where monomers with two attractive bonding sites 
polymerize into polydisperse rods). Moreover, one 
often encounters phrases in the literature, such 
as "This theory [Onsager's theory] shows that re- 
pulsive interactions [excluding volume interactions] 
alone can lead to long-range orientational nematic 
order, disproving the notion that attractive inter- 
actions are a prerequisite" [17] , which could be am- 
biguous with respect to the role that attractive lat- 
eral interactions between the rods should play in 
reinforcing (or not) the nematic order. 

In this context, it is of interest and of value to 
inquire how the existence of lateral interactions be- 
tween the rods influences the phase transition oc- 
curring in the system. The objective of this paper 
is to provide a thorough analysis in this direction. 
For this purpose, an exhaustive study of the phase 
transition occurring in a system of attractive rigid 
rods deposited on square lattices was performed. 
The results revealed that (i) the orientational or- 



der survives in the presence of attractive lateral in- 
teractions; (ii) the critical density shifts to higher 
values as the magnitude of the lateral interactions 
is increased; and (Hi) the continuous transition be- 
comes first order for interaction strength w > w c 
(in absolute values). 

The outline of the paper is as follows. In Sec. 
II we describe the lattice-gas model and the sim- 
ulation scheme. In Sec. Ill we present the MC 
results. Finally, the general conclusions are given 
in Sec. IV. 

II. Lattice-gas model and Monte 
Carlo simulation scheme 

We address the general case of adsorbates assumed 
to be linear rigid particles containing k identical 
units (fe-mers), with each one occupying a lattice 
site. Small adsorbates would correspond to the 
monomer limit (fc = 1). The distance between fc- 
mer units is assumed to be equal to the lattice con- 
stant; hence exactly k sites are occupied by a k-mer 
when adsorbed (see Fig. [1]). The surface is repre- 
sented as an array of M = L x L adsorptive sites 
in a square lattice arrangement, where L denotes 
the linear size of the array. In order to describe 
the system of N fc-mers adsorbed on M sites at a 
given temperature T, let us introduce the occupa- 
tion variable Cj which can take the values Cj = 
if the corresponding site is empty and Cj = 1 if 
the site is occupied. On the other hand, molecules 
adsorb or desorb as one unit, neglecting any pos- 
sible dissociation. Under these considerations, the 
Hamiltonian of the system is given by 

H = w CjCj - N(k - l)w + e ^ (1) 

where w is the nearest-neighbor (NN) interaction 
constant which is assumed to be attractive (nega- 
tive), represents pairs of NN sites and e is 
the energy of adsorption of one given surface site. 
The term N(k — l)w is subtracted in eq. ([I]) since 
the summation over all the pairs of NN sites over- 
estimates the total energy by including N(k — 1) 
bonds belonging to the N adsorbed fc-mers. Be- 
cause the surface was assumed to be homogeneous, 
the interaction energy between the adsorbed dimcr 
and the atoms of the substrate e Q was neglected for 
the sake of simplicity. 
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Figure 1: Linear tetramers adsorbed on square lat- 
tices. Full and empty circles represent tetramer 
units and empty sites, respectively. 

In order to characterize the phase transition, we 
use the order parameter defined in Ref. [llj . which 
in this case can be written as 

6= lnh - Uvl (2) 
rih + n v 

where rih(n v ) is the number of rods aligned along 
the horizontal (vertical) direction. When the sys- 
tem is disordered (6 < 9 C ), all orientations are 
equivalents and 6 is zero. As the density is in- 
creased above Q Cl the fc-mers align along one direc- 
tion and S is different from zero. Thus, 5 appears 
as a proper order parameter to elucidate the phase 
transition. 

The problem has been studied by grand canon- 
ical MC simulations using a typical adsorption- 
desorption algorithm. The procedure is as fol- 
lows. Once the value of the chemical potential 
fi is set, a linear fc-uple of nearest-neighbor sites 
is chosen at random and an attempt is made to 
change its occupancy state with probability W — 
min{l,exp(-AiJ/fc B T)}, where AH = H f -H l is 
the difference between the Hamiltonians of the final 
and initial states and fcs is the Boltzmann constant. 
In addition, displacement (diffusional relaxation) of 
adparticles to nearest-neighbor positions, by either 
jumps along the fc-mer axis or reptation by rotation 
around the fc-mer end, must be allowed in order to 
reach equilibrium in a reasonable time. A MC step 
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Figure 2: Adsorption isotherms (coverage versus 
chemical potential) for k = 10, L = 100 different 
i/j/igT's as indicated. Inset: Adsorption phase di- 
agram of attractive 10-mers on square lattices. 

(MCs) is achieved when M fc-uples of sites have 
been tested to change its occupancy state. Typi- 
cally, the equilibrium state can be well reproduced 
after discarding the first r' = 10 7 MCs. Then, the 
next r = 2 x 10 7 MCs are used to compute averages. 

In our MC simulations, we varied the chemical 
potential n and monitored the density 9 and the 
order parameter 5, which can be calculated as sim- 
ple averages. The reduced fourth-order cumulant 
Ul introduced by Binder [18] was calculated as: 

^ = 1 -^y2' ( 3 ) 

where (• • • ) means the average over the MC simu- 
lation runs. All calculations were carried out using 
the BACO parallel cluster (composed by 60 PCs 
each with a 3.0 GHz Pcntium-4 processor and 90 
PCs each with a 2.4 GHz Core 2 Quad processor) 
located at Instituto de Ffsica Aplicada, Universi- 
dad Nacional de San Luis-CONICET, San Luis, 
Argentina. 

III. Results 

The calculations were developed for linear 10-mers 
(k = 10). With this value of k and for non- 
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interacting rods, it is expected the existence of a 
nematic phase at intermediate densities [10] . The 
surface was represented as an array of adsorptive 
sites in a square lattice arrangement with conven- 
tional periodic boundary conditions. The effect of 
finite size was investigated by examining lattices 
with L = 50,100,150,200. 

In order to understand the basic phenomenology, 
we consider, in the first place, the behavior of the 
adsorption isotherms in presence of attractive lat- 
eral interactions between the A:-mers. 

Fig. [2]shows typical adsorption isotherms (cover- 
age versus ^t/ksT) for linear 10-mers with different 
values of the lateral interaction (the solid circles 
represent the Langmuir case, w/kBT = 0). 

The isotherms shift to lower values of chemical 
potential, and their slopes increase as the ratio 
w/ksT increases (in absolute value). For inter- 
action strength above the critical value (w > w c , 
in absolute values) the system undergoes a first- 
order phase transition, which is observed in the 
clear discontinuity in the adsorption isothermfl In 
the case studied, this critical value is approximately 
w c /k B T w -0.80 (or k B T c /w w -1.25). The be- 
havior of the adsorption isotherms also allows us to 
calculate the phase diagram of the adsorbed mono- 
layer in "temperature-coverage" coordinates. In 
fact, once obtained the real value of the chemical 
potential (or critical chemical potential /x c ) in the 
two-phase region, the corresponding critical densi- 
ties can be easily calculated. By repeating this pro- 
cedure for different temperatures ranging between 
and T c , the coexistence curve can be built [20] , 
A typical phase diagram, obtained in this case for 
attractive 10-mers, is shown in the inset of Fig. O 

On the basis of the study in Fig. [2 our next 
objective is to obtain evidence for the existence of 
nematic order in the range —0.80 < w/ksT < of 
attractive interactions. For this purpose, the be- 
havior of the order parameter 5 as a function of 
coverage was analyzed for k — 10, L — 100 and dif- 
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Figure 3: Surface coverage dependence of the ne- 
matic order parameter for k — 10, L = 100 different 
u>//cbT"s as indicated. 



ferent values of the lateral interaction. The results 
are shown in Fig. [3J revealing that (i) the orien- 
tational order survives in the presence of attrac- 
tive lateral interactions and (n) the critical density 
shifts to higher values as the magnitude of the lat- 
eral interactions is increased. 

In order to corroborate the results obtained in 
the last figure, we now study the dependence of 
8 C on w/ksT. In the case of the standard the- 
ory of FSS [TSl [TO], when the phase transition is 
temperature driven, the technique allows for vari- 
ous efficient routes to estimate T c from MC data. 
One of these methods, which will be used in this 
case, is from the temperature dependence of Ul{T), 
which is independent of the system size for T = T c . 
In other words, T c is found from the intersection 
of the curve Ul(T) for different values of L, since 
Ul{T c ) =const. In our study, we modified the con- 
ventional FSS analysis by replacing temperature by 
density [TT]. Under this condition, the critical den- 
sity has been estimated from the plots of the re- 
duced four-order cumulants Ul{0) plotted versus 
6 for several lattice sizes. As an example, Fig. [4] 
shows the results for w/ksT = —0.125. In this 
case, the value obtained was 9 C — 0.542(2). In 
the inset, the data are plotted over a wider range 
of temperatures, exhibiting the typical behavior of 
the cumulants in the presence of a continuous phase 
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Figure 4: Curves of U L {6) vs for k = 10, 
w/ksT = —0.125 and square lattices of different 
sizes. From their intersections one obtained 9 C . In 
the inset, the data are plotted over a wider range 
of densities. 



Figure 5: Temperature-coverage phase diagram 
corresponding to attractive fc-mers with k 
I ii. The inset in the upper- left (lower-right) cor- 
ner shows a typical configuration in the nematic 
(isotropic) phase. 



transition. 

The procedure of Fig. [4] was repeated for 
—0.80 < w/ksT < 0, showing that the values of 
8 C increase linearly with the magnitude of the lat- 
eral couplings (see solid squares in Fig. [5]). The 
critical line (dotted line in the figure) was obtained 
from the linear fit of the numerical data. As it is 
possible to observe, the range of coverage at which 
the transition occurs diminishes as w/ksT is in- 
creased (in absolute value). This finding indicates 
that the presence of attractive lateral interactions 
between the rods does not favor the formation of 
nematic order in the adlayer. The phenomenon 
can be understood from the behavior of the sec- 
ond virial coefficient, which will initially decrease 
on introducing attractive w. This decrease implies 
that the isotherms shift to lower values of chemical 
potential, and consequently, the critical point shifts 
to higher densities. 

We did not assume any particular universality 
class for the transitions analyzed here in order to 
calculate their critical densities, since the analy- 
sis relied on the order parameter cumulant's prop- 
erties. However, the fixed value of the cumu- 
lants, U* — 0.617(9), is consistent with the ex- 
tremely precise transfer matrix calculation of U* = 



0.6106901(5) [21] for the 2D Ising model. This find- 
ing may be taken as an indication that the phase 
transition belongs to the 2D Ising universality class. 

With respect to the behavior of the system for 
w/ksT < —0.80, the adsorbed layer "jumps" from 
a low-coverage phase to a high-coverage phase. 
This effect, which has been discussed in Fig. [2j 
is represented in Fig. [5] by the dashed coexistence 
line. The low-coverage phase is an isotropic state, 
similar to that observed for w/ksT > —0.80 and 
low density (see inset in the lower-right corner of 
Fig. [5]). On the other hand, the high-coverage 
phase is also an isotropic state, but characterized by 
the presence of local orientational order (domains 
of parallel fc-mers). A typical configuration in this 
regime is shown in Fig. [5]). 

Finally, it is worth pointing out that: (1) the be- 
havior of the order parameter in Fig. [3] clearly indi- 
cates that the transition from the low-density dis- 
ordered phase to the intermediate-density ordered 
phase is an isotropic to nematic phase transition 
(when all the words have the usual meaning). In 
this case, the transition under study belongs to the 
2D Ising universality class. It can also be thought 
of as an unmixing or liquid-gas transition |15j . For 
this reason we have called gas and liquid to the 
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Figure 6: Typical configuration of the adlayer in 
the high-coverage phase and w/IcbT < —0.80. 

phases reported in Fig. and (2) even though it 
has not been rigorously proved yet, a second phase 
transition for non-interacting rods at high densi- 
ties has been theoretically predicted [10] and nu- 
merically confirmed [13]. This result has not been 
confirmed for the case of attractive rods. An ex- 
haustive study on this subject will be the object of 
future work. 

IV. Conclusions 

We have addressed the critical properties of attrac- 
tive rigid rods on square lattices with two allowed 
orientations, and shown the dependence of the crit- 
ical density on the magnitude of the lateral interac- 
tions w/ksT. The results were obtained by using 
MC simulations and FSS theory. 

Several conclusions can be drawn from the 
present work. On the one hand, we found that 
even though the presence of attractive lateral in- 
teractions between the rods does not favor the for- 
mation of nematic order in the adlayer, the orien- 
tational order survives in a range that goes from 
w/ksT = up to Wc/ksT « —0.80 (w c //csT rep- 
resents the critical value at which occurs a typical 
transition of condensation in the adlayer). In this 
region of w/ksT, the critical density increases lin- 
early with the magnitude of the lateral couplings. 
On the other hand, the evaluation of the fixed point 
value of the cumulants U* — 0.617(9) indicates 



that, as in the case of non-interacting rods, the ob- 
served phase transition belongs to the universality 
class of the two-dimensional Ising model. 

With respect to the behavior of the system for 
w/ksT < —0.80, the continuous transition be- 
comes first order. Thus, the adsorbed layer jumps 
from a low-coverage phase, similar to that observed 
for w/ksT > —0.80 and low density, to an isotropic 
phase at high coverage, characterized by the pres- 
ence of local orientational order (domains of paral- 
lel fc-mers) 

Future efforts will be directed to (1) extend the 
study to repulsive lateral interaction between the 
A:-mers; (2) obtain the whole phase diagram in 
the space (temperature-coverage- rod's size); (3) de- 
velop an exhaustive study on critical exponents and 
universality and (4) characterize the second phase 
transition from a nematic order to a non-nematic 
state occurring at high density. 
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